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Abstract 

The concept of "self-organized criticality" (SOC) has been introduced by Bak, 
Tang, and Wiesenfeld (1987) to describe the statistics of avalanches on the surface 
of a sandpile with a critical slope, which produces a scale-free powerlaw size distri- 
bution of avalanches. In the meantime, SOC behavior has been identified in many 
nonlinear dissipative systems that are driven to a critical state. On a most general 
level, SOC is the statistics of coherent nonlinear processes, in contrast to the Pois- 
son statistics of incoherent random processes. The SOC concept has been applied 
to laboratory experiments (of rice or sand piles), to human activities (population 
growth, language, economy, traffic jams, wars), to biophysics, geophysics (earth- 
quakes, landslides, forest fires), magnetospheric physics, solar physics (flares), 
stellar physics (flares, cataclysmic variables, accretion disks, black holes, pulsar 
glitches, gamma ray bursts), and to galactic physics and cosmology. 
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1 Introduction 

Physical processes in our universe can be subdivided into incoherent and coherent pro- 
cesses. Examples of incoherent processes are fractional Brownian motion, diffusion, 
collisional plasma processes, plasma heating, which all are governed by random pro- 
cesses and thus can be described by binomial and Poissonian statistics, which exhibit 
exponential-like or Gaussian-like distributions. In contrast, coherent processes over- 
come the threshold of the random noise background and grow in a multiplicative way, 
such as avalanches, chain reactions, or catastrophes, which all exhibit powerlaw-like 
size distributions. The dychotomy of incoherent and coherent processes is also re- 
flected in the two types of linear and nonlinear processes. Linear systems are character- 
ized by a proportionality between the input and output. Small disturbances cause small 



1 




Figure 1: The output or dissipated energy in a linear system grows linearly with time, for a 
constant input rate (top left), while the output is highly unpredictable and not correlated with the 
input rate in a nonlinear dissipative system (bottom left). A practical example of a linear system 
is a hydroelectric plant, where the produced electric energy is proportional to the water input, as 
depicted with the water-storage dam at Yaotsu, Gifu, Japan (top right). A classical example of a 
dissipative nonlinear system is a snow avalanche, as shown in the large wet-snow avalanche at 
Deadman Canyon in the Sierra Nevada range (bottom right). 

effects, and large disturbances are required to produce a large effect. A hydropower 
plant, for instance, produces electric energy that is proportional to the water inflow 
rate and gravitational or kinetic energy of the water that is feeding the turbine (Fig.l, 
top). Nonlinear processes, in contrast, are governed by a multiplicative or amplification 
factor, that leads to outputs of unpredictable magnitude, unlike linear systems. Small 
disturbances can cause both small or large avalanches. Classical examples are snow 
avalanches (Fig. 1, bottom), landslides, floodings, forest fires, or earthquakes. 

A conceptual model to understand the basic nature of nonlinear processes was in- 
troduced by Bak, Tang, and Wiesenfeld (1987) in terms of sandpile avalanches. We can 
build a sandpile by dripping sand grains randomly with a more or less steady rate, until 
it forms a cone (Fig. 2, right). Once the sandpile reaches a critical slope (typically at an 



2 




SUB-CRITICAL STATE : CRITICAL STATE : 

NO AVALANCHES AVALANCHES 



Figure 2: A static equilibrium produces no avalanche events (bottom left panel), such as the flat 
sand beach in northern California (top left panel), while randomly dripping sand onto a sandpile 
produces a state of self-organized criticality where avalanches occur (bottom right panel), such 
as with the conveyer belt of the Indian River Enterprises (top right panel). 



angle of repose of « 34° — 37°), avalanches occur in an unpredictable way, even when 
the input of dripped sand grains is steady. Bak called this driven state " self-organized 
criticality" (SOC), which has become a paramount characteristic of many nonlinear 
systems. A sand beach, in contrast, would be in a sub-critical state without avalanches 
(Fig. 2, left). The key features of nonlinear systems are powerlaw-like size distribu- 
tions, which indicate scale-free parameter distributions, such as the size, time scale, or 
energy of an avalanche. The occurrence frequency distribution N(x) of a parameter x 
is defined as a powerlaw function when 

N(x) oc x~ a . (1) 

The mathematical nature of such powerlaw distributions can easily be understood in 
terms of a coherent, exponential growth process as described in the following section. 
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Figure 3: Time evolution of energy release rate Wit) for 3 different ratios of growth times 
to saturation times, rajts = (0.5, 1.0,2.0) (left) and the corresponding powerlaw distribu- 
tions of the peak energy release rate P. Note that the event set with the shortest growth time 
(rafts = 0.5) reaches the highest energies and thus produces the flattest powerlaw slope 

(a = 1 + T G /t s = 1.5). 



2 Analytical Model of SOC Avalanches 

Avalanches occurring in the state of self-organized criticality represent local instabil- 
ities that grow explosively for some time interval. The released energy grows in a 
nonlinear way above some energy threshold, which can be parameterized by some 
nonlinear function, for instance by an exponential growth function. We define the time 
evolution of the energy release rate W{t) of a nonlinear process that starts at a threshold 
energy of Wq by 



W(t) = W exp 



< t < t , 



(2) 



where tq represents the exponential growth time. The process grows exponentially 
until it saturates at time t = t with a saturation energy Ws, 

W s = W(t = t) = W cxp (j-^J . (3) 

We define a peak energy release rate P that represents the maximum energy release rate 
Ws, after subtraction of the threshold energy Wo, that corresponds to the steady-state 
energy level before the nonlinear growth phase, 



P = W s - Wo = Wo 



(4) 



In the following, we will refer to the peak energy release rate P also briefly as "peak 
energy". For the saturation times r, which we also call "rise times", we assume a 
random probability distribution, approximated by an exponential function N(t) with 
e-folding time constant ts, 



N(T)dT = — exp 
ts 



ts) 



(5) 
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This probability distribution is normalized to the total number of N events, 

N(r)dT = N . (6) 



/ 
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In order to derive the probability distribution N(P) of peak energy release rates P, we 
have to substitute the variable of the peak energy, P, into the function of the rise time 

t(P), 

dr 

N(P)dP = N(r)dT = N[t(P)} 



dP 

This requires the inversion of the evolution function P(t) (Eq. 4), 



dP . (7) 



and the calculation of its derivative dr/dP, which is 

dr tq ( P ' 



r(P) = t g In ( £ + 1 ) , (8) 



dP Wo \W ' 7 ' ^ 

Inserting the probability distribution of saturation times N(t) (Eq. 5), the inverted 
evolution function t(P) (Eq. 8) and its time derivative {dr/dP) from Eq. (9) into the 
frequency distribution N(P) in Eq. (7) yields then, 

which is an exact powerlaw distribution for large peak energies (P 3> Wo) with a 
powerlaw slope a of 

*=(i + f). mi 

The same mathematical result was derived in Rosner and Vaiana (1978), but the ex- 
ponential growth was attributed to energy storage therein, rather than to the nonlinear 
growth phase of the instability here. The powerlaw slope thus depends on the ratio 
of the growth time to the e-folding saturation time, which is essentially the average 
number of growth times. Examples of time series with avalanches of different growth 
times (rafts — 0.5, 1.0, 2.0) are shown in Fig. 3, along with the corresponding pow- 
erlaw distributions of peak energies P. Note that the fastest growing events produce 
the flattest powerlaw distribution of peak energies. 

This simple model explains the powerlaw distribution of peak energies in nonlinear 
processes that have an exponential-like growth, which is the key characteristic of co- 
herent processes. In incoherent processes, the growth time would be much larger than 
the average saturation time (tq ~> t$), in which case the time evolution (Eq. 3) can 
be linearized, leading to P oc t (Eq. 4), and the distribution of peak energies would be 
identical to that of the saturation times t, which is a Poissonian (exponential) distribu- 
tion (Eq. 5) with a characteristic scale P s=y Wq. So, incoherent processes would have 
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Figure 4: Time evolution of total released energy e(t) (top panels), the energy release rate 
W(t) — de(t)/dt (middle panels), and binary representation of avalanche growth rate (bottom 
panels), for both the exponential (left panels) and the logistic growth model (right panels). An 
exponential curve (right top) and a Gaussian curve (right middle) are drawn (with dotted lines) 
onto the logistic curves for comparison (Aschwanden et al. 1998). 



the same exponential distribution of peak fluxes P as their random durations r, 

N^{P)dP=£ o ev(-£r)dP. (12) 

Therefore, we can use the functional shape of the occurrence frequency distribution 
of event sizes as a diagnostic of linear or nonlinear processes: exponential distribution 
functions indicate an incoherent random process, while powerlaw distribution func- 
tions indicate a coherent growth process. 

The exponential-growth model represents just the mathematically simplest function 
to produce a powerlaw distribution of peak energies. It corresponds to a multiplicative 
process where the number of avalanche elements doubles every time step until it satu- 
rates (Fig. 4, left). The saturation phase, however, abruptly ends with a discontinuity. 
A more natural model is the logistic equation (discovered by Pierre Francois Verhulst 
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Figure 5: Solar flares observed in EUV with the TRACE spacecraft in 171 A: The flare of 2001 
Apr 15 exhibits an erupting filament in the foreground and a rising postflare arcade behind near 
the limb (left panel), while the 2000 Nov 9 flare displays the 3D geometry of the double-ribbon 
postflare arcade (right panel) [courtesy of NASA/TRACE]. 



in 1845), which is a first-order differential equation that approaches the saturation limit 
asymptotically, 



dE(t) _ E(t) 



dt 



TG 



1 



E(t) 



(13) 



as shown in Fig. 4 (right). The time derivative, which mimics the energy release rate, 
has a maximum when the total energy has the steepest slope. The avalanche grows 
first exponentially and dies out gradually after the peak time of the fastest growth. It 
can be shown mathematically that the frequency distribution of the peak energies is 
also a powerlaw (Aschwanden et al. 1998), and thus this more physical description of a 
logistic avalanche with a continuous time evolution shares the same characteristic as an 
exponentially-growing avalanche, and thus both models can be used to describe SOC 
processes. 



3 Solar Observations of SOC Phenomena 

Solar flares are probably the best-studied datasets regarding SOC statistics in astro- 
physics. Solar flares are catastrophic events in the solar corona, most likely caused by a 
magnetic instability that triggers a magnetic reconnection process, producing emission 
in almost all wavelengths, such as in gamma-rays, hard X-rays, soft X-rays, extreme 
ultraviolet (EUV), Ha emission, radio wavelengths, and sometimes even in white light. 
Two EUV images of solar flares are shown in Fig. 5, which display the intricated topol- 
ogy of the magnetic field produced at the beginning (Fig. 5, left) and end phase of the 
flare (Fig. 5, right). Since the emission mechanisms are all different in each wave- 
length, such as nonthermal bremsstrahlung (in hard X-rays), thermal bremsstrahlung 
(in soft X-rays and EUV), gyrosynchrotron emission (in microwaves), plasma emis- 
sion (in metric and decimetric waves), etc., we expect that the calculation of energies 
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Figure 6: Occurrence frequency distributions of hard X-ray peak count rates P(cts/s) ob- 
served with HXRBS/SMM (1980-1989), BATSE (1991-1993), and RHESSI (2002-2007), with 
powerlaw fits. Note that BATSE/CGRO has larger detector areas, and thus records higher count 
rates. 

contained in each event strongly depends on the emission mechanism, and thus on the 
wavelength. It is therefore advisable to investigate the statistics of SOC events in each 
wavelength domain separately. The most unambiguous SOC parameters to report are 
the peak flux P, the total flux or fluence E, defined as the time-integrated flux over 
the entire event, and the total time duration T of the event. Conversions of fluxes and 
fluences into energy release rates and total energies require physical models. 

Statistics of the peak count rates P of over 20,000 solar flares observed in hard 
X-rays during three different space missions (SMM 1980-1989; CGRO 1991-1993; 
RHESSI 2002-2010) and during three different solar cycles have been gathered and 
their occurrence frequency distributions (or log N - log S distributions) are shown in 
Fig. 6. All three datasets show a near-perfect powerlaw distribution function with a 
mean slope of ap = 1.75 ± 0.05. According to our model this value corresponds to 
a ratio of rafts = 0.75 (Eq. 11) or a mean saturation time ts/ra = 1.33 growth 
times. The distributions shown in Fig. 6 span together over 5 orders of magnitude, 
which means that the largest flare had a maximum saturation time of In (10 5 ) = 11.5 
growth times, which corresponds to an extremely coherent process that could never be 
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Figure 7: Composite flare frequency distribution in a normalized scale in units of 1CP 50 flares 
per time unit (s -1 ), area unit (cm -2 ), and energy unit (erg -1 ). The energy is defined in terms 
of thermal energy E t h = 3n e kBT e V for EUV and SXR, and in terms of non-thermal energy 
in > 25 keV (Crosby et al. 1993) or > 8 keV electrons (Lin et al. 2001). The slope of -1.8 is 
extended over the entire energy domain of 10 24 — 10 32 erg. The offset between the two HXR 
datasets is attributed to different lower energy cutoffs as well as different levels of flare activity 
during the observed time intervals [adapted from Aschwanden et al. 2000 and Lin et al. 2001]. 

explained by random probability. 

Similar powerlaw-like occurrence frequency distributions have also been found for 
other flare parameters, such as the fluence, the flare duration, or flare energy (calculated 
from the thermal energy or energy in nonthermal hard X-ray producing electrons). 
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Figure 8: Cumulative occurrence frequency distributions 7V(> E) of stellar flares (vertically 
shifted for each star, in order of stellar distance) are shown as a function of the flare energy 
E, observed from 12 different stars by Audard et al. (2000). For comparison we show also 
the statistics of total hard X-ray emission of 18,888 solar flares observed with RHESSI, which 
covers an energy range of E « 10 27 — 10 31 ' 7 ergs, while stellar flares have an energy range of 
E ~ 10 30 ' 7 — 10 35 - er g S Note that most distributions show a (cumulative) powerlaw part with 
a slope of a « 0.8 (diagonal lines) at the low end and an exponential high-energy drop-off. 

Based on this fact of scale-free powerlaw distributions, Lu and Hamilton (1991) were 
the first to apply the concept of self-organized criticality (SOC) to solar flares. With a 
cellular automaton model that mimics sandpile avalanches, they were able to simulate 
the size frequency distribution with a powerlaw slope of ap = 1.67 ± 0.02, which is 
close to the observed value of ap = 1.75 ± 0.05 shown in Fig. 6. 

The statistics of peak energies in solar flares has been extended to other wave- 
lengths, from hard X-rays to soft X-rays and EUV and the powerlaw-like size distribu- 
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Figure 9: Frequency distribution of giant-pulse flux densities measured from the Crab pulsar, 
observed during 15-27 May 1991 with the Green Bank 43-m telescope at 1330, 800, and 812.5 
MHz. The tail can be represented by a powerlaw distribution Nf oc F~ a with a = 3.46 ± 0.04 
for fluxes F > 200 Jy (Lundgen et al. 1995). 

tion has been found to extend over 8 orders of magnitude (Fig. 7). The events observed 
in soft X-rays have been called "active region transients", having about 3-5 orders of 
magnitude less energy than the largest solar flare. The events observed in EUV have 
energies about 6-9 orders of magnitude less than the largest flare, and thus are called 
"nanoflares". The self-similarity of flare parameters observed over this wide wave- 
length range, which implies also a large temperature range of T w (1 — 35) x 10 6 
K, constitutes a powerful argument that the solar corona is in a state of self-organized 
critic ality. 

4 Astrophysical Observations of SOC Phenomena 

Size distributions have also been sampled for stellar flares, although with much smaller 
statistics, typically < 15 flare events per star. A comparison of size distributions from 
12 flare stars (type F to M) observed by Audard et al. (2000) is shown in Fig. 8. 
Because of the smallness of the samples, the cumulative frequency distributions are 
shown, which have slope that is flatter by one compared with the differential fre- 
quency distribution (ct c p m = otp — 1). The comparison in Fig. 8 demonstrates two 
important aspects: (1) stellar flares are up to two orders of magnitude more energetic 
(E « 10 32 - 10 34 ) than solar flares (see RHESSI comparison at bottom of Fig. 8), 
although the small sample implies common and relatively weak events, and (2) the 
drop-off at the upper bound of the cumulative frequency distribution yields a steeper 
slope than the powerlaw slope obtained from more complete frequency distributions. 
Correcting for the latter effect, stellar flares seem to have a similar powerlaw distribu- 
tion as solar flares, and thus may be goverend by the same SOC process. 
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Other size distributions with powerlaw shape have been observed in astrophysical 
data forpulsar glitches (e.g., Lundgren et al. 1995; shown in Fig. 9), for soft gamma-ray 
repeaters (e.g., Gogus et al. 1999), accretion disk objects around black hole candidates, 
such as Cygnus X-l (e.g., Negoro et al. 1995), or forblazars (e.g., Ciprini et al. 2003). 
These astrophysical objects display size distributions with different powerlaw values, 
and thus can be explained by different physical processes, but the distribution type of a 
powerlaw function constitutes a powerful argument that the underlying system is driven 
into the critical state of self-organized criticality, which is common to earthquakes, 
sunquakes (occurring during large flares), and neutron star quakes (pulsar glitches). 
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